Bethe[L_, M_] := (
unknowns = Table[u[k], {k,M}];
eqns=Table[
(((((u[k]+I/2)/(u[k]-I/2))))^L)==-(∏_(j=1)^M)((u[k]-u[j]+I)/(u[k]-u[j]-I)),
{k,M}
];
sols=Select[
Solve[eqns, unknowns],
Length[Union[Last /@ #]]==M&
];
Es=2(∑_(k=1)^M)(1/(((u[k])^2)+1/4)) /. sols;
{eqns, sols, Es // N}
)

Bethe[4,0]

Bethe[4,1]

Bethe[4,2]

Bethe[4,3]

Bethe[5,0]

Bethe[5,1]